function q = Elnpmu(mix)

% E[ln p(mu)]

[K D] = size(mix.centres);
Dlog2pi = D*log(2*pi);

m0 = mix.varprior.m0;
inv_T0 = inv(mix.varprior.T0);
logdet_T0 = logdet_chol(mix.varprior.T0);

m = mix.varposterior.m;
T = mix.varposterior.T;

q = 0;
for k = 1:K
    q = q - 0.5*Dlog2pi - 0.5*logdet_T0 - 0.5*(m(k, :)-m0)*inv_T0*(m(k, :)-m0)' - 0.5*trace(T(:, :, k)*inv_T0);
end
